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ABSTRACT 

The pulsar wind model is updated by considering the effect of particle density and 
pulsar death. It can describe both the short term and long term rotational evolution of 
pulsars consistently. It is applied to model the rotational evolution of the Crab pulsar. 
The pulsar is spun down by a combination of magnetic dipole radiation and particle 
wind. The parameters of the Crab pulsar, including magnetic field, inclination angle, 
and particle density are calculated. The primary particle density in acceleration region 
is about 10 3 times the Goldreich-Julian charge density. The lower braking index be¬ 
tween glitches is due to a larger outflowing particle density. This may be glitch induced 
magnetospheric activities in normal pulsars. Evolution of braking index and the Crab 
pulsar in P — P diagram are calculated. The Crab pulsar will evolve from magnetic 
dipole radiation dominated case towards particle wind dominated case. Considering 
the effect of pulsar “death”, the Crab pulsar (and other normal pulsars) will not evolve 
to the cluster of magnetars but downwards to the death valley. Different acceleration 
models are also considered. Applications to other sources are also discussed, including 
pulsars with braking index measured, and the magnetar population. 

Key words: pulsars: general - pulsars: individal (PSR B0531+21) - stars: magnetar 
- stars: neutron - wind 


1 INTRODUCTION 

The Crab pulsar (PSR B0531+21) is a young radio pulsar 
with a spin frequency v = 30.2 Hz and frequency derivative 
z> = —3.86 x 10~ lo Hz/s (Lyne et al. 1993). Its character¬ 
istic magnetic field is about 7.6 x 10 12 G at the magnetic 
polcfl The Crab pulsar has been monitored continuously 
for decades years by various telescopes since it was discov¬ 
ered in 1968 (Lyne et al. 1993; Lyne et al. 2015; Wang et 
al. 2012). Its braking index is n = 2.51 ± 0.01 (Lyne et al. 
1993). Different braking index values 2.45, 2.57 and 2.3 are 
also reported between glitches (Wang et al. 2012; Lyne et 
al. 2015). Observational details are given in Table [T] 

Long-term observations have found that almost all pul¬ 
sars (except accreting X-ray pulsars in binary systems) are 
spinning down. The spin-down behavior can be described by 
a power law (Lyne et al. 1993): 

z> = —kv n , (1) 

here, k is usually taken as a constant and n is the brak- 


1 Assuming all the rotational energy loss is due to magnetic 
dipole radiation in vaccum, B c = 6.4 X 10 19 \JPP G 


ing index. The braking index and second braking index are 
defined respectively (Livingstone et al. 2005): 


( 2 ) 

(3) 


where v is the pulsar spin frequency, z>, v and v are the 
first, second and third frequency derivative, respectively. 
The spinning down of pulsars is usually assumed to be 
braked by the magnetic dipole radiation (Shapiro & Teukol- 
sky 1983). In the magnetic dipole radiation model, a neutron 
star rotates uniformly in vacuo at a frequency v and pos¬ 
sesses a magnetic moment /x. The corresponding slowdown 
rate is: 


8ttV 3 . 2 

v = —n# v sm Q ’ 


(4) 


where /x = 1/2 BR 3 is the magnetic dipole moment (B is the 
polar magnetic field and R is the radius of neutron star), 
/ = 10 45 g • cm 2 is the moment of inertia, c is the speed of 
light, and a is the angle between the rotational axis and the 
magnetic axis (i.e., the inclination angle). The spin down 
behavior of pulsar in this model can be described as v oc 
zx 3 . The braking index is exactly three if /x, /, and a are 
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constant. To date, only eight pulsars, including the Crab 
pulsar, have measured the meaningful braking indices for 
they are young and own relatively larger v (Espinoza et al. 
2011; Lyne et al. 2015). Their braking indices are all smaller 
than three. It means that there are other physical processes 
needed to slow down the pulsar. 

Many mechanisms have been proposed in order to ex¬ 
plain the braking index observations, e.g., the pulsar wind 
model (Xu & Qiao 2001; Wu et al. 2003; Contopoulos & 
Spitkovsky 2006; Yue et al. 2007), a changing magnetic field 
strength (Lin et al. 2004; Chen & Li 2006; Espinoza et al. 
2011), a changing inclination angle (Lyne et al. 2013), and 
additional torques due to accretion (Liu et al. 2014). How¬ 
ever, these models are more or less not consistent with obser¬ 
vations or can not simulate the long-term evolution of pul¬ 
sars. And these models can not explain the different braking 
indices detected between glitches (Wang et al. 2012; Lyne 
et al. 2015). Furthermore, the effect of pulsar death (death 
line, Ruderman & Sutherland 1975; or death valley, Chen 
& Ruderman 1993; Zhang et al. 2000) should be considered 
in modeling the long term rotational evolution of the pulsar 
(Contopoulos & Spitkovsky 2006). 

The pulsar wind model considers both the pulsar spin- 
down and the particle acceleration in the magnetosphere 
(Xu & Qiao 2001; Wu et al. 2003). In this paper, an updated 
pulsar wind model is built based on previous researches (Xu 
& Qiao 2001; Yue et al. 2007; Contopoulos & Spitkovsky 
2006; Li et al. 2012). It includes: (1) Both magnetic dipole 
radiation and particle outflow, and their dependence on the 
inclination angle; (2) The particle outflow depends on the 
specific acceleration model; (3) The primary particle density 
may be much larger than the Goldreich- Julian charge den¬ 
sity; (4) The effect of pulsar death is considered in modeling 
the rotational evolution of the pulsar. This model can calcu¬ 
late both the short term and long term evolution of pulsars. 
It is applied to the Crab pulsar which has the most detailed 
timing observations. Possible applications to other sources 
are also illustrated. 

The pulsar wind model and model calculations of the 
Crab pulsar are presented in section 2. Discussions and con¬ 
clusions are given in section 3 and section 4, respectively. 


2 SPIN-DOWN OF THE CRAB PULSAR IN 
THE PULSAR WIND MODEL 


2.1 Description of the pulsar wind model 


In the pulsar wind model, the rotational energy is consumed 
by magnetic dipole radiation and particle acceleration (Xu & 
Qiao 2001). The magnetic dipole radiation is related to the 
perpendicular component of the magnetic moment (p,± = 
^tsina), the power of magnetic dipole radiation is (Shapiro 
& Teukolsky 1983): 


E d 


2/r 2 H 4 

3c 3 


• 2 

sm a, 


(5) 


where H = 2 nu is the angular speed of the pulsar. The effect 
of parallel component of magnetic moment (/in = /rcosa) is 
responsible for particles acceleration (Ruderman & Suther¬ 
land 1975). Acceleration gaps are formed above the polar 
gap, primary particles are generated and accelerated in these 
gaps. The energy taken away by these particles dependents 


Table 1. Observations of the Crab pulsar. 


Pulsar parameter Values 


Epoch(MJD) 

t'(Hz) 

i>(10- 10 Hz/s) 
i>(10 _20 Hz/s 2 ) 
v (10 — 30 Hz/s 3 ) 

Braking index n 
Second braking index m 
Inclination angle ct(°) 
Age 


40000.0 

30.225437“ 

-3.86228“ 

1.2426“ 

-0.64“ 

2.51 ± 0.01“ 

10.15“ 

(45 ~ 70) 6 

915 yr (from 1054 to 1969) 


Notes: (a): Observed spin frequency and its derivatives are 
parameters in the Taylor expansion 
v = uq + vo(t — to) + \i'o(t - to) 2 + | i/o (t — to) 3 for 
1969 ~ 1975 (Lyne et al. 1993). For there are rare glitches 
occurred in this interval. Later observations (Lyne et al. 2015) 
confirm previous results. 

(b):The range of inclination angle is given by modeling the pulse 
profile of the Crab pulsar (Dyks et al. 2003; Harding et al. 2008; 
Watter et al. 2009; Du et al. 2012; Lyne et al. 2013). 


on the acceleration potential drop. The corresponding rota¬ 
tional energy loss rate is (Xu & Qiao 2001): 

E p = 27rrpCp e A0, (6) 

where r p = R(Rflfc) 1 ^ 2 is polar gap radius, p e = «Pgj is the 
primary particle density (pGJ = HR/(27rc) is the Goldreich- 
Julian charge density, Goldreich & Julian 1969 ), A <f> is the 
corresponding acceleration potential of the acceleration re¬ 
gion, and k is a coefficient related to primary particle density 
which can be constrained by observation^. The maximum 
acceleration potential for a rotating dipole is A<b = pH 2 /c 2 
(Ruderman & Sutherland 1975). Then equation <j6j> can be 
rewritten as (considering that the particle acceleration is 
mainly related to the parallel component of magnetic mo¬ 
ment): 


E P 


2p 2 H 4 

3c 3 


o A0 2 
COS 

A<f> 


a. 


(7) 


The magnetic dipole radiation and the outflow of particle 
wind may contribute independently. Then the total rota¬ 
tional energy loss rate is: 


E 


■ • 2 u 2 n 4 2 Ad> 2 , 

Ed + Ep = (sm a + 3k —— cos a) 

3c 3 A<f> 

2 p 2 H 4 


( 8 ) 


where p = sin 2 a + 3kA<(>/A4> cos 2 a. The first and second 
items of expression p are respectively for magnetic dipole 
and particle wind. If the acceleration potential A <j> = 0, there 
are no particles accelerated in the gap, the pulsar is just 
braked down by the magnetic dipole radiation. The accel¬ 
eration potential is model dependent (Xu & Qiao 2001; Wu 
et al. 2003). A particle density of k times Goldreich-Julian 


2 Primary particle density in the acceleration gap is: p e = 
|e|[?r_|_ + n_], but Goldreich-Julian “charge” density is: pgj = 
|e|[n_|_ — n_] (Yue et al. 2007). It is reasonable to infer that: 

K > 1. 
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charge density is considered. Expressions of r) are given in 
Table[2]for various acceleration modeL0 The rotational evo¬ 
lution of the pulsar can be written as: 


-/nn = 


2 n 2 n 4 

oTs r l’ 


(9) 


According equation © and ©. the braking index and sec¬ 
ond braking index in the wind braking model are: 


Q. dr) 

n = 3+ - — , 
77 ail 


m = 


K , d V , dr). 


( 10 ) 

( 11 ) 


Equation @, disj and m are all functions of magnetic 
held B, inclination angle a and particle density k (fl and its 
derivatives are observational inputs, see Table [l]). 

Taking these expressions of r) (Table [2J back to equa¬ 
tion ©. it can be seen that: as the pulsar spinning down (fl 
decreases), the effect of particle wind becomes stronger. If 
there are particles flowing out, it means that there are still 
particle acceleration in magnetosphere (the pulsar is still ac¬ 
tive). However, as pulsar spinning down, the potential drop 
may not sustain the need to accelerate particles (Ruderman 
& Sutherland 1975). Besides, observations do not support 
pulsars to evolve unceasingly in these pulsar wind models 
(Young et al. 1999; Contopoulos & Spitkovsky 2006). There¬ 
fore, “pulsar death” must be considered in modeling the long 
term rotational evolution of the pulsar. As the spin angu¬ 
lar speed decreasing to the death value, the radio emission 
tends to stop. And when pulsar is death, the pulsar is braked 
only by magnetic dipole radiation. Here, the VG(CR) model 
is employed to show the effect of pulsar death by introduc¬ 
ing a piecewise function deduction (details are shown in the 
appendix A). 


CR 


IVG = 


sin 2 a + 4.96 X 10 2 x(l - n ~ 15 / 7 cos 2 a, 

if > 0 ^ eat h ( 12 ) 

sin 2 a , 

if Q < 


The death period is defined as fldeath = 27r/Pdeath (Con¬ 
topoulos & Spitkovsky 2006; Tong & Xu 2012): 

Pdeath = 2.8( I ^) 1/2 ( I ^)- 1/2 S. (13) 

V sap can be viewed as the maximum acceleration poten¬ 
tial drop in the open field line regions and 14ap = 10 13 V 
is chosen in the following calculations (Contopoulos & 
Spitkovsky 2006). It is introduced according to Contopoulos 
& Spitkovsky (2006) (see also later works by Li et al. 2012). 
In this way, the effect of pulsar death can be incorporated 
in the rotational energy loss rate. The effect of pulsar death 
is discussed in more detail in section 2.5. For young pulsars 
(fi 3> fldeath), the effect of pulsar death is negligible. But for 
the long-term evolution of pulsars, especially when the pe¬ 
riod approaches the death period, the effect of pulsar death 
is significant. 


3 As shown in Table[2] the models of VG(CR) and SCLF(I) are 
very similar with each other. Crab parameters calculated in these 
two models are also similar (shown later). 


Table 2. The Expressions of r) for nine particle acceleration mod¬ 
els. 


No. Acceleration model 77 


1 VG (CR) 

2 VG (ICS) 

3 SCLF (II,CR) 

4 SCLF (II, ICS) 

5 SCLF (I) 

6 OG 

7 CAP 

8 NTVG (CR) 

9 NTVG (ICS) 


sin 2 ol + 4.96 X 10 2 kB^ 8 ^ 7 Cl 15 / 7 cos 2 ol 

sin 2 a + 1.02 x 10 5 kB~ 2 2 ^ 7 Q ~ 13 / 7 cos 2 ol 

sin 2 ol + 38kB~ 2 Cl~ 7 / 4 cos 2 ol 

sin 2 ol + 2.3/cS 12 22//13 f2 —8 / 13 cos 2 a 

sin 2 ol + 9.8 X 10 2 K,B 12 8/7 n~ 15 / 7 cos 2 ol 

sin 2 a + 2.25 X 10 5 kB” 12/7 0 _26/7 cos 2 a: 

sin 2 ol + 54 k,B~ 2 Cl~ 2 cos 2 a 

sin 2 ol + 13.7/3° ^kB~^ f2 —1-76 cos 2 a 

• 2 1 — 1 D —10 — 1-88 2 

sin ct + 69.67 k ^12 ^ cos a 


Notes: Particle density p e = repQJ i s taken in all these models. B 12 is the 
magnetic field strength in units of lO ^ 2 G. 

(a) : The models 1 ~ 5 are based on the pulsar wind model of Xu 8z Qiao 
(2001). The VG and SCLF are respectively the acceleration models: 
vacuum gap and space charge limit flow. In the vacuum gap (Ruderman & 
Sutherland 1975), curvature radiation (CR) and inverse Compton 
scattering (ICS) are considered (Zhang et al. 2000). In the SCLF case, 
regimes II and I are defined by field saturated or not (Arons &c 
Scharlemann 1979; Harding & Muslimov 1998). Three models are 
introduced: CR for SCLF model in regime II, ICS for SCLF model in 
regime II, and the SCLF model in regime I. 

(b) : The OG means the outer gap, is self-sustaining and limited by the 
electron-positron pair produced by collisions between high-energy photons 
(Zhang &z Cheng 1997). The modification by Wu et al. (2003) is adopted. 

(c) : The CAP model is a phenomenological model of constant potential 
A<£ = 3 X 10 12 V (Yue et al. 2007). 

(d) : The NTVG shorts for near threshold vacuum gap. In this model, the 
electron-positron pair produced at or near the kinematic threshold 

hu = 2mc 2 / sin 9 because of superstrong surface magnetic field (Gil &; 
Melikidze 2002). /3 = 52 is taken in the CR case, and 7 = 14 is taken in 
the ICS case (Abrahams & Shapiro 1991; Wu et al. 2003). 


2.2 Understanding the braking index of the Crab 
pulsar 

Parameters of the Crab pulsar are calculated in this section. 
The VG(CR) model is taken as an example to show the cal¬ 
culation process. For pulsars with the observed v, v, v and 
v, their observational braking index n and second braking 
index m can be get by equation m and ©. Pulsar param¬ 
eters such as magnetic field, inclination angle, and particle 
density can be calculated by these three equations ©, G 0 
and (HU. The primary particle density may be 10 3 ~ 10 4 
times of Goldreich-Julian charge density (Yue et al 2007). 
According the observational range of inclination angle and 
the characteristic magnetic field (which can be viewed as 
order of magnitude estimation of the true magnetic field), 
the range of particle density k can be limited by equation 
© and m- In our calculation, k = 10 3 is adopteejf], the 
inclination angle a and magnetic field B are about 55° and 
8.1 x 10 12 G which are consistent with observational con¬ 
straints. The second braking index calculated by equation 
dm is 10.96, which is roughly consistent with the observed 
value 10.15 (Lyne et al. 1993), considering the observational 
uncertaintities (Lyne et al. 1993, 2015). 

The magnetic field and inclination angle are assumed 
to be constant in the pulsar wind model. The particle den¬ 
sity decreases in long-term evolution of pulsars but remains 
unchanged in short duration at early time. Figure |T| shows 
the braking index of the Crab pulsar as a function of spin 


4 We should note that the k is related to the primary particles 
in the acceleration gap but not the total out-flow particles. 
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2.3 An increasing particle density results in a 
lower braking index during glitches 



Figure 1. The braking index of the Crab pulsar as function of 
spin period in the VG(CR) model. The dashed line is observa¬ 
tional braking index of 2.51. The dotted line is braking index in 
the transition point (n = 2). 


Table 3. Parameters of the Crab pulsar calculated in various 
acceleration models. 


Models 

K 

10 3 

O' 

0 

&12 

Po 

ms 

Pt 

ms 

Pd 

s 

T t 

10 3 yr 

T d 

10 5 yr 

VG(CR) 

1 

55 

8.1 

18.8 

55 

2.52 

2.77 

0.78 

VG(ICS) 

0.1 

60 

7.5 

18.9 

63 

2.58 

3.50 

1.67 

SCLF(II CR) 

2 

58 

7.5 

18.9 

68 

2.59 

3.93 

2.24 

SCLF(II ICS) 

0.6 

44 

5.0 

19.2 

- 

3.16 

- 

30.1 

SCLF(I) 

0.6 

57 

7.8 

18.8 

55 

2.45 

2.77 

0.78 

OG 

20 

59 

8.2 

18.5 

42 

2.47 

1.63 

0.07 

CAP 

3 

52 

8.3 

18.9 

58 

2.46 

3.08 

1.14 

NTVG(CR) 

3 

57 

7.7 

19.5 

67 

2.56 

4.10 

2.52 

NTVG(ICS) 

30 

54 

7.5 

19.4 

62 

2.58 

3.58 

1.81 


Notes: 

(a) : k is the coefficient of primary particle density p e = kpqj. 

(b) :a and B 12 are respectively the inclination angle and the magnetic 
field. 

(c) : Pq is the initial period. 

(d) :Pt and Tt are the transition period and corresponding age. At this 
transition point the braking index is 2. Because the minimum braking 
index of model SCLF(ICS) is 2.4, there is no such a transition point, 
which can be seen in the evolution in P — P diagram (see Figure HI. 

(e) : Pj and Tj are the period and age when pulsar stops radio emission. 
The death period is calculated by equation ( i 13 ft . 


period. As the pulsar period increases, the braking index 
decreases from 3 to 1. It illustrates that the Crab pulsar is 
initially braked by magnetic dipole (n 3) and then by 
the particle wind (n —> 1). We take n = 2 as the transition 
point (see Figure [5] at n = 2, the slop of the curve is 0) and 
the corresponding spin period is about 55 ms. Since Table 
|T] has shown the timing observations of the Crab pulsar in 
1969 (its age is 915 years old by then), its initial period is 
Po fei 18.8 ms in AD 1054 by integrating equation ©. The 
calculations of all acceleration models are shown in Table [3] 
As shown in this table, the particle densities are different 
in different models, but are all more than 100 times of pgj- 
The inclination angles is within the observational range. The 
initial periods are all about 19 ms. 


A braking index 2.3 has been monitored after the removal of 
the effects of glitches during an epoch when there are many 
glitches occured (Lyne et al. 2015). It is different from the 
long term underlying braking index 2.51 of the Crab pulsar. 
Previously, Wang et al. (2012) have measured two different 
values of braking index n = 2.45 and n = 2.57 between 
glitches. They attribute it to the effect of varying particle 
wind strength (Wang et al. 2012). It is generally accepted 
that glitch is caused by the inner effect but may lead to 
some effect in the outer magnetospher^l, for example: the 
changing of primary particle density during this epoch. In 
the pulsar wind model, it can be denoted by the changing 
coefficient of particle density k = K,{t). Then equation (1101) 
should be modified: 

Q. dri k dri ft k , 

n = 3 H--|- -L -—. 14 

i] di l rj dK Q k 

Equation m can be rewritten as: 


n drj 

n = 3 + --777 
77 di l 


k dr/ r c 
r) dnrT 


(15) 


where r c = — is the characteristic age, t k = can be 
viewed as the typical time scale of the changing particle den¬ 
sity. When the magnetosphere is in equilibrium (i.e., when 
the glitch activities of the pulsar is not very active), the 
time scale of particle density variation (t k ) may be very 
large (i.e., larger than the characteristic age r c ). The third 
term in equation m does not contribute. But, if they are 
comparable with each other, the braking index changes sub¬ 
stantially. I 11 other words, when the outflow particle density 
increases and the effect of particle wind becomes stronger, 
means t k > 0 or k > 0, the braking index will be smaller 
than 2.51. But when the outflow particle density decreases 
(t k < 0 or k < 0), the braking index will be larger than 
2.51. As the out-flow particle density tending to a certain 
value (larger or smaller than previous value), the braking 
index tends to its underlying value (2.51). Generally, the 
increased (or decreased) component of the outflow particle 
density (may be 0.1%) is so much small that it can be ig¬ 
nored. 

Figure [2] shows the braking index as function of r K for 
the Crab pulsar in the VG(CR) model. I 11 order to bet¬ 
ter understand the observations of n < 2.51 (Lyne et al. 
2015), we consider the r K > 0 only. As we can see in this 
figure, the braking index is insensitive to r K when it is larger 
than 10 4 yr. But, when it is comparable with the character¬ 
istic age (about 10 3 yr), the braking index decreases sharply. 
When braking index is 2.3, the changing rate of the parti¬ 
cle density is k « 1.2 x 10 -8 s _1 . In this interval (from MJD 
51000 to MJD 53000, about 5 years), the particle density has 
changed by 2 pgj • The change of particle density will result 
in changes of radiation energy loss rate and period deriva¬ 
tive. Their changes respectively are: SE « 2 x 10 35 erg/s and 
8P « 1.88 x 10 -16 s/s. 


5 Glitch induced magnetospheric activities are very common in 
the case of magnetars (Dib &; Kaspi 2014). 
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Figure 2. The braking index of the Crab pulsar as function of r K 
in the VG(CR) model. The dashed line is the underlying braking 
index 2.51. The dotted line is the smaller braking index 2.3 mea¬ 
sured during glitch activities (Lyne et al. 2015, Figure 7 there). 


The pulsar wind model can be tested with observations 
(Jodrell Bank monthly ephemeris of the Crab pulsai0). The 
effect of glitches is not taken into consideration in the pul¬ 
sar wind model. The changes of period and period deriva¬ 
tive caused by glitches should be added (Lyne et al. 2015). 
Figure [3] is the model calculation compared with observa¬ 
tions (without modelling the transient variation caused by 
glitches). When the amplitude of glitches is relatively small, 
the change of particle density is not obvious, then the model 
calculations (the blue points) fit well with observations (the 
black points). However, when the amplitude of glitches is 
large, the effect of glitches must be added (the red points). 
Moreover, the effect of particle density change should be 
taken into consideration. A factor SP = 1.88 x ICC 16 s/s 
(which is calculated above) is added to the red points (from 
MJD 51804.75 afterwards), shown as the green points. The 
green points can match the general trend of the observations. 


2.4 Long term evolution of the Crab pulsar 

The pulsar period and its fist derivative evolution with time 
in the VG(CR) model are shown in the first and second 
panel of Figure [I] The curves evolve slowly in its early age 
and then rise sharply, indicating that the Crab pulsar is 
mainly braked by magnetic dipole radiation firstly and then 
mainly by particle wind. The transition period is 55 ms and 
age of the Crab pulsar will be 2771 calculated in the VG(CR) 
model (see Table O. The braking index evolution with time 
in the VG(CR) model is shown in the third panel of Figure 
U Braking index decreases from 3 to 6/7 -different mini¬ 
mum braking indices for different acceleration models. The 
bottom panel of Figure 0 shows the second braking index 
evolution with time. The curve gradually decreases from 15 
to 30/49. Figure [5] shows the evolution of the Crab pul¬ 
sar in P — P diagram from birth to its 30000 years old. 
These points are respectively taken when the Crab pulsar 
is: 1, 915, 2771, 9000 and 30000 years old. Period derivative 


6 http://www.jb.man.ac.uk/pulsar/crab.html, data from 1982 
February 15 (MJD 45015) to 2014 October 15 (MJD 56945) are 
used here. 



Figure 3. Evolution of the Crab pulsar in P — P diagram in 
VG(CR) model compared with observations. The black points 
are timing results from Jodrell Bank monthly ephemeris. The 
blue points are the rotational evolution of Crab pulsar in the 
VG(CR) model. The red ones are obtained by adding the effect of 
glitches (Lyne et al. 2015, Table 3 there). And the green ones are 
the summation of red points and a factor 5P ~ 1.88 X 10 1 s/s 
(from MJD 51804.75 afterwards), which may be caused by the 
change of out-flow particle density. At first, the black, blue, red, 
and green points are coincide with each other. (Notes: The fitting 
does not model the transient process caused by glitches.) 

evolves with period by a function of P oc P 2 ~ n (Espinoza 
et al. 2011). In the log-log plot, the evolution curve evolves 
with a slop (2 — n). As the braking index falling from 3 to 1, 
the pulsar evolve from magnetic dipole radiation dominated 
case (the left points) to the particle wind dominated case 
(the right points) and the bottom is the transition point 
(n = 2). Clearly, the curve evolves more sharply after the 
transition point. The asymptotic behaviors are discussed in 
the appendix B. 


2.5 The effect of pulsar death 

As shown is the Figure [5] the line evolves up right with a 
slope about 1 after the transition point just like PSR J1734- 
3333 whose braking index is 0.9±0.2 (Espinoza et al. 2011). 
If the pulsar continually evolve in this trend, it would not 
dead and might end itself as a magnetar. It will be hard 
to explain the observations for: (1) The number of magne- 
tars are smaller than radio pulsars and most of the pulsars 
occupy the central region in the P-P diagram (see Figure 
O; (2) Only a small portion of magnetars are radio emitters 
(Olausen & Kaspi 2014); (3) The properties of the Crab 
pulsar, PSR J1734-3333 and other high magnetic field pul¬ 
sars are significantly different from magnetars (Ng & Kaspi 
2011). As pulsars slow-down, the density of outflow particles 
may decrease and the effect of particle wind will recede (for 
the SCLF cases, the number of particles (or k) may close 
to constant, but the particle density still decreases because 
a reducing Goldreich-Julian charge density). Therefore, the 
effect of pulsar death must be included in modeling the long¬ 
term rotational evolution of pulsars. 

The essential condition of radio emission for a pulsar is 
its pair production (Sturock 1971; Ruderman & Sutherland 
1975). The so called pulsar “death” means the stopping of 
pulsar emission (Ruderman & Sutherland 1975; Chen & Ru¬ 
derman 1993; Zhang et al. 2000; Contopoulos & Spitkovsky 
2006). The death defined by Ruderman & Sutherland (1975) 
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Figure 4. Rotational evolution of the Crab pulsar from birth to 
30 thousands years in the VG(CR) model. The first and second 
panels are respectively the period and period derivative evolution 
with time. The points are observations of period and its derivative 
(Lyne et al. 1993). The third and fourth panels are respectively 
the braking index and second braking index evolution with time. 
The dashed lines are observed values and the solid lines are re¬ 
spectively the minimum braking index and second braking index 
in the VG(CR) model. 
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Figure 5. Evolution of the Crab pulsar in the P-P diagram. Five 
points are chosen when the pulsar is: 1, 915, 2771, 9000 and 30000 
years old. 



is that: a pulsar dead when the maximum potential drop A® 
available from the pulsar is smaller than the acceleration 
potential needed to accelerate particles. Because this death 
line is model dependent, Chen & Ruderman (1993) proposed 
death “valley” by modifying the boundary conditions. And 
this death “valley” was updated by considering different 
particle acceleration and photon radiation processes (Zhang 
et al. 2000). However these works are just to separate pul¬ 
sars with radio emission from those without. Contopoulos & 
Spitkovsky (2006) included the effect of pulsar death when 
modeling the rotational evolution of pulsars. However, the 
braking indices there are always larger than three. Consid¬ 
ering particle acceleration and pulsar death simultaneously 
may give a comprehensive interpretation for both short and 
long-term rotational evolution of pulsars. 

The visual explanation for pulsar death in the pulsar 
wind model is that the particles in the magnetosphere are 
exhausted and there are no primary particles generated. 
For the physical understanding, we can refer to the equa¬ 
tion m . expressions A<3> oc f l 2 and Atf> (i.e., A()>vgcr = 
2.8 x 10 13 R6 /7 -Br 2 1/7 O- 1/7 V, Xu & Qiao 2001), with the 
pulsar spinning down, the maximum potential of the pulsar 
(A<F) drops, when the maximum potential can not meet the 
need to accelerate particles (A tj> > A4>), the pulsar death. 
Besides, the acceleration potential is weakly dependent on 
the fi, that is why we can make calculations under the con¬ 
stant acceleration potential (the CAP model) (Yue et al. 
2007). Figure [6] shows the long-term evolutions of the Crab 
pulsar in the VG(CR) model. The dashed curve is the line 
shown in Figurc[5]and the solid is the one with pulsar death 
(according to equation (1121) 1. These two lines evolve simi¬ 
lar at the early age, and divide obviously after the transi¬ 
tion point. The second line falls down afterwards and moves 
down-right with a slop about —1 finally. The initial period 
and transition period can be calculated according the second 
line, which are same with the values given by the first line 
within precision. Indicating that the effect of pulsar death 
can be ignored in the early time. As mentioned before, the 
particle density decreases as pulsar spin-down. The gap be¬ 
tween these two lines is caused by the reducing of particle 
density. When the spin period of Crab approaches to the 
death period 2.52 s, the radio emission tends to stop and 
the slop of the second line is very large. Crab evolution in 
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Figure 6. Rotational evolution of the Crab pulsar in the VG(CR) 
model. The dashed curve is the evolution without considering pul¬ 
sar death. The solid line is the one considering pulsar death. The 
P — P diagram has listed all known radio pulsars, magnetars, 
and millisecond pulsars (updated from Figure 1 in Tong & Wang 
2014). The dot-dashed line represents the fiducial death line (Ru- 
derman &: Sutherland 1975). Pulsars with meaningful braking 
indices (Espinoza et al. 2011; Lyne et al. 2015) are defined by 
large square. The arrows represent their evolution directions and 
each arrow indicate its motion in the next 10000 yr (For the Vela 
pulsar, 20000 yr is used; 5000 yr for PSR J1846-0258 and PSR 
J1119-6127). 

various acceleration models are given in Figure [7] The bot¬ 
tom line is the evolution curve in SCLF(ICS) model. There 
is no transition point, indicating that the particle wind ef¬ 
fect of this model is very weak. The top line is the evolution 
of the Crab pulsar in the OG model. This model is for high- 
energy radiation, the death of this model may differ from 
radio emission models. The results here for the OG model 
are just crude approximations to show the effect of pulsar 
death. Evolution of the Crab pulsar in other models are sim¬ 
ilar with each other. 


3 DISCUSSIONS 

In the magnetic dipole radiation model, a pulsar is assumed 
as a magnetic dipole rotating in vacuum. It is just a low- 
order approximation of the true magnetosphere, and there 
must be particles acceleration and generation process so that 
the radio emission can be received. The actual magneto¬ 
sphere surrounding the pulsar must be filled with plasma 
and corotate with pulsar because of magnetic freeze. Accord¬ 
ing the special relativity, the corotation magnetosphere can 
not infinitely extend but only exist within the light cylin¬ 
der. The magnetosphere is divided into two parts by the 
light cylinder: the closed magnetic field lines and the open 
magnetic field lines (including polar gap and outer gap). In 
the open field line region, particles are accelerated, gener¬ 
ate radiation and flow out (thus taken away some of the 
rotational energy of the central neutron star). In the pulsar 


Figure 7. Rotational evolution of the Crab pulsar in all accel¬ 
eration models. Because the VG(CR) model and SCLF(I) model 
are very similar, their lines are coincident. Only the line in the 
VG(CR) model is plotted. See also Figure[6]for explanations. 

wind model, the pulsar is braked down by magnetic dipole 
radiation and particle wind. But in the electric current loss 
model, the electric current flow generates a spin-down torque 
(Harding et al. 1999; Beskin et al. 2007). Based on the same 
magnetospheric physics, the pulsar wind model and the cur¬ 
rent loss case are just different point of views, will produce 
the same results (e.g. Section 3.2 in Tong et al. 2013). And 
they can be used to check each other. It is meaningful to take 
advantage of the existing information and explain more ob¬ 
servations. The particle wind worked on pulsar spin-down 
was first proven by observations of intermittent pulsar PSR 
B1931+24 whose radio emission switches between “on” and 
“of” state (Kramer et al. 2006). Corresponding calculations 
for intermittent pulsars in pulsar wind model were made by 
Li et al. (2014). 

As one of the best observed pulsar, observations of the 
Crab pulsar are relatively comprehensive (see Table|Tj& tim¬ 
ing results from Jodrell Bank website). It will be the best 
target to test radiation models. The pulsar parameters can 
be calculated by equations ©, (E3 and GD when model de¬ 
pendent acceleration potential (or r/) is given (see Table [2|- 
The Goldreich-Julian charge density is taken as the primary 
particle density in these models of Xu & Qiao (2001) and Wu 
et al (2003). However, Yue et al. (2007) constrained the pri¬ 
mary particle density by observations of several pulsars and 
predicted that it should be 10 3 ~ 10 4 times of Goldreich- 
Julian charge density. In our calculations, the primary par¬ 
ticle density is at least 100 times of Goldreich-Julian charge 
density. In the VG(CR) model, a particle density 10 3 |pgj| is 
chosen. When applying the pulsar wind model to intermit¬ 
tent pulsars, a Goldreich-Julian charge density is assumed 
(Li et al. 2014). As the pulsar spin-down, its out-flow particle 
density decreases. In the P — P diagram, these intermittent 
pulsars are on the right of those young pulsars (see Figure 
©. This may explain the different particle density in the 
Crab pulsar and intermittent pulsars. 
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Other studies also found the need for a high particle 
density in the pulsar magnetosphere, e.g., when modeling 
the eclipse of the double pulsar system, giant pulses of the 
Crab pulsar, and magnetar X-ray spectra (Lyutikov 2008 
and references therein), optical/ultraviolet excess of X-ray 
dim isolated neutron stars (Tong et al. 2010, 2011). Mag- 
netohydrodynamical simulations of the Crab nebula also re¬ 
quire a much higher electron flow from the magnetosphere 
(Buhler & Blandford 2014 and references therein). Besides, 
a super Goldreich-Julian current has been proposed in a se¬ 
ries theoretical works about the pairs creations in the accel¬ 
eration gap (Timokhin & Arons, 2013; Zanotti et al. 2012; 
Hirotani 2006). Combined with the results in this paper, it 
is possible that there are high density plasmas (e.g., 10 3 -10 4 
times the Goldreich-Julian density) in both open and closed 
field line regions in the pulsar magnetosphere. 


3.1 Applications to other sources 

Here the updated pulsar wind model is employed to com¬ 
pare with the observations of the Crab pulsar. It can also 
be applied to other sources. For eight pulsars with meaning¬ 
ful braking index measured (Espinoza et al. 2011; Lyne et 
al. 2015), their magnetospheric parameters can also be cal¬ 
culated, including magnetic field, inclination angle, particle 
density etc. But the higher second frequency derivatives of 
some pulsars (Hobbs et al. 2010) may be attributed to the 
fluctuation of magnetosphere (Contopoulos 2007). During 
this process, more information will be helpful to constrain 
the model parameters, e.g., inclination angle, second braking 
index. Their evolution on the P — P diagram can be calcu¬ 
lated, which will be similar to Figure [6] The pulsar PSR 
J1734—3333 has braking index n = 0.9 ± 0.2 (Espinoza et 
al. 2011). Its rotational energy loss rate will be dominated 
by the particle wind. Also from Figure [6] the braking indices 
of young pulsars will naturally be divided into two groups. 
Braking indices of pulsars in the first group are close to 
three, their evolution directions are down to the right in the 
P — P diagram (P oc P 2 ~ n ~ P _1 ). Braking indices in the 
second group are close to one and their evolution directions 
are up to the right (P oc P 2 ~ n ~ P 1 )- Generally speaking, 
sources in the second group will be older than the ones in 
the first group. At present, there are some hints for the evo¬ 
lution of pulsar braking index (Espinoza 2013). Future more 
observations will deepen our understanding on this aspect. 

Glitch is a phenomenon of the pulsar suddenly spin¬ 
ning up and then slowing down by a larger rate. In the pul¬ 
sar wind model, the lower braking index 2.3 during glitches 
(of the Crab pulsar) is due to the larger out-flow particle 
density. The increasing particle density is a possible perfor¬ 
mance associated with glitches. Glitches induced magneto- 
spheric activities are very common in magnetars (Kaspi et 
al. 2003; Dib & Kaspi 2014). For the high magnetic field pul¬ 
sar PSR J1846—0258, it also shows a smaller braking index 
after glitch (Livingstone et al. 2011). The cause of a smaller 
braking index may be similar to the Crab pulsar case. Due to 
a larger particle outflow, the pulsar should experience some 
net spindown after glitch. This has already been observed 
(Livingstone et al. 2010). Some primary calculations for PSR 
J1846—0258 has been done previously (Tong 2014, based on 
a wind braking model designed for magnetars). The smaller 


braking index of PSR J1846—0258 after glitch is consistent 
with the results here. 

From Figure [6] and Figure [7] the pulsar will go up-right 
in the P — P diagram in the wind braking dominated case. 
Therefore, for pulsars in the up right corner of the P — P 
diagram (e.g., high magnetic field pulsars and magnetars), 
their true dipole magnetic field may be much lower than the 
characteristic magnetic field. For magnetar SGR 1806—20, 
its characteristic magnetic field is about 5 x 10 15 G (Tong 
et al. 2013). This dipole magnetic field is too high to under¬ 
stand comfortably (Vigano et al. 2013). Using the physical 
braking mechanism in this paper, the dipole magnetic field 
of SGR 1806—20 will be much lower. In the late stage of a 
pulsar, it will bent and go downward in the P — P diagram. 
The Crab pulsar and other pulsars will not evolve into the 
cluster of magnetars. Furthermore, for magnetars, they also 
will go downward in the P — P diagram when they are aged 
enough. This may corresponds to the case of low magnetic 
field magnetars (Tong & Xu 2012). 


4 CONCLUSIONS 

An updated pulsar wind model which considered the effect of 
particle density and pulsar death is developed. It is employed 
to calculate parameters and simulate the evolution of the 
Crab pulsar. The braking index n < 3 is the combined effect 
of magnetic dipole radiation and particle wind. The lower 
braking indices measured between glitches are caused by the 
increasing particle density. By adding the glitch parameters 
and a change of period derivative caused by the change of 
particle density, the theoretical evolution curve is well fitted 
with observations (Figure [3J. This may be viewed as glitch 
induced magnetospheric activities in normal pulsars. Giving 
the model dependent acceleration potential drop, the mag¬ 
netic field, inclination angle and particle density of the Crab 
pulsar can be calculated. The evolution of braking index and 
the Crab pulsar in P — P diagram can be obtained. In the 
P—P diagram, the Crab pulsar will evolve towards the death 
valley, not to the cluster of magnetars (Figure [6]). Different 
acceleration models are also considered. The possible appli¬ 
cation of the present model to other sources (pulsars with 
braking index measured, and the magnetar population) is 
also mentioned. 


ACKNOWLEDGMENTS 

The authors would like to thank the Referee Barsukov D.P. 
for helpful comments and R.X.Xu, J.B.Wang for discussions. 
H.Tong is supported Xinjiang Bairen project, West Light 
Foundation of CAS (LHXZ201201), Qing Cu Hui of CAS, 
and 973 Program (2015CB857100). 


APPENDIX A: CONSIDERATION OF PULSAR 
“DEATH” 

The treatment of Contopoulos & Spitkovsky (2006) is em¬ 
ployed in the following. Assuming that the open field line 
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regions rotate at a angular velocity n open but not the angu¬ 
lar velocity of the neutron stafl 


APPENDIX B: ASYMPTOTIC BEHAVIROS IN 
FIGURE 4 AND 5 
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The particle density in this acceleration gap is k times the 
Goldreich-Julian charge density: 
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The energy taken away by the out-flowing particles is: 


E p = 2nrlcp e A(j>. 


(A3) 


The perpendicular component and parallel component of the 
magnetic moment are respectively related to the magnetic 
dipole radiation and particle acceleration. The rotation en¬ 
ergy loss rate can be written as: 
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with the expression of r / 
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For a constant acceleration potential Ac/) = 3 X 10 12 V, cor¬ 
respondingly, 


rj = sin 2 a + 54k(1 — 
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For the physical acceleration models, r/ is model dependent. 
Its expression in the VG(CR) case is calculated in the fol¬ 
lowing. The potential drop of open field lines is (Ruderman 
& Sutherland 1975): 


A<t>= QopenB h 2 , (A7) 

c 

where f2 open is the angular frequency of the open held 
lines (Contopoulos & Spitkovsky 2006), h is the thick¬ 
ness of the gap, h = 1.1 x 10 4 p^ 7 cm and 

p = 2.3 x 10 8 cm is the curvature radius of the 

open held lines (pe = p/10 6 cm) (Ruderman & Sutherland 
1975). The expression of r/ can be rewritten as: 

p = sin 2 a + 4.96 x 10 2 k 

x (l - ^g^) 6/7 R- 19/7 R- 8/7 «“ 15/7 cos 2 a (A8) 

The exponent of (1 — n<1 ^ ath ) is the minimum braking index 
in the VG(CR) model, n = 6/7 (Li et al. 2014). In the wind 
braking model, the minimum braking index is always around 
one (Li et al. 2014). Therefore, it is taken as one in equation 
(1 1211 (for all the acceleration models). During the numerical 
calculations, the neutron star radius is taken as 10 6 cm (i.e., 

R& = !)• 


7 In fact, particles in the acceleration gap rotate with a angular 
speed between the spinning speed of pulsar surface and the open 
field lines (Ruderman & Sutherland 1975). We use the angular 
speed of the open field lines as the boundary condition to describe 
the pulsar death because the acceleration potential is insensitive 
to it in this pulsar wind model. 


In the pulsar wind model, the pulsar is braked down by the 
combination of magnetic dipole radiation and particle wind. 
The effect of particle wind is more and more important with 
pulsar spinning down. As shown in Figured the period and 
period derivative evolve up sharply with age. If the pulsar 
braking is dominated by magnetic dipole radiation, equation 
CD can be written as: 

v = -kv 3 . (Bl) 

the evolution of spin frequency is v = (2kt)~ 1 ^ 2 (assuming 
initial spin frequency is very large). Correspondingly, the 
period and period derivative evolution with time are respec¬ 
tively: P = (2 kt) 1 ^ 2 , and P = k(2kt)~ 4 ^ 2 . If the pulsar is 
braked down mainly by particle wind, the braking index is 
about 1. The spin-down power law is : 

v = —kv. (B2) 

the evolution of spin frequency is v oc exp (— kt). Corre¬ 
spondingly, the period and period derivative evolution with 
time are respectively: P oc exp(fcf), and P oc exp(fcf). 
Clearly, the period and period derivative evolve faster by 
the exponential form. 
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